Influence of gas compression on flame acceleration in the early stage of burning in 
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Abstract 

The mechanism of finger fiame acceleration at the early stage of burning in tubes was studied experimentally by Clanet and 
Searby [Combust. Flame 105: 225 (1996)] for slow propane-air flames, and elucidated analytically and computationally by 
Bychkov et al. [Combust. Flame 150: 263 (2007)] in the limit of incompressible flow. We have now analytically, experimentally 
and computationally studied the finger fiame acceleration for fast burning fiames, when the gas compressibility assumes an 
important role. Specifically, we have first developed a theory through small Mach number expansion up to the first-order terms, 
demonstrating that gas compression reduces the acceleration rate and the maximum flame tip velocity, and thereby moderates 
the finger fiame acceleration noticeably. This is an important quantitative correction to previous theoretical analysis. We 
have also conducted experiments for hydrogen- oxygen mixtures with considerable initial values of the Mach number, showing 
finger flame acceleration with the acceleration rate much smaller than those obtained previously for hydrocarbon flames. 
Furthermore, we have performed numerical simulations for a wide range of initial laminar flame velocities, with the results 
substantiating the experiments. It is shown that the theory is in good quantitative agreement with numerical simulations for 
small gas compression (small initial flame velocities). Similar to previous works, the numerical simulation shows that finger 
fiame acceleration is followed by the formation of the "tulip" flame, which indicates termination of the early acceleration 
process. 

Keywords: premixed flames, flame acceleration, compressibility, finger flames, tulip flames 
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sound speed 

heat capacity at constant pressure 

heat capacity at constant volume 

activation energy 

Lewis number 

flame thickness 

flame length in a 2D configuration 

initial flame propagation Mach number 

pressure 

Prandtl number 

energy diffusion vector 

energy release from the reaction 

tube radius (channel half-width) 

universal gas constant 

Schmidt number 

unstretched laminar burning velocity 

time 

temperature 

velocity 

dimensionless axial velocity 

dimensionless radial velocity 

radial coordinate 

mass fraction of fuel 

axial coordinate 

auxiliary constant 

adiabatic index 

total energy per unit volume 

stress tensor 

dimensionless radial coordinate 



e 

p 

a 

T 
TR 



gas expansion ratio 
instantaneous expansion factor 
dynamical viscosity 
dimensionless axial coordinate 
density 

scaled acceleration rate 
scaled time (r = S'^i/i?) 
factor of time dimension 



Subscripts and other designations 






initial 


pi 


planar geometry 


axi 


axisymmetric geometry 


b 


burnt gas 


C 


modified 


f 


flame, flame skirt 


s 


sidewall 


sj 


just ahead of the flame skirt 


s,b,f 


in the burnt gas at the flame skirt 


sph 


instant/locus of shape change 


tip 


flame tip 


u 


unburnt mixture 


w 


wave 


wall 


flame touches side walls of tube 



scaled value 



1. INTRODUCTION 

The spontaneous acceleration of a flamefront prop- 
agating from the closed end of a tube is a key ele- 
ment in deflagration-to-detonation transition (DDT) 01- 
[S], with two main mechanisms identified as possible 
causes for the flame acceleration, namely the mechanisms 
of Shelkin [6j and Bychkov et al. [7j for flame propaga- 
tion in smooth and obstructed tubes respectively. Ac- 
cording to the Shelkin mechanism, flames accelerate in 
smooth tubes due to the non-slip boundary conditions at 
the walls. Quantitative theory of the process has been de- 
veloped and validated by extensive numerical simulations 
in Refs. [HIS]- For the Bychkov mechanism, the delayed 
burning between the obstacles in obstructed tubes pro- 
duces a powerful jet flow driving an extremely fast flame 
acceleration in the unobstructed, central portion of the 
tube [71 [10]. 

In addition to these two major mechanisms, another 
possible scenario of flame acceleration in tubes has been 
demonstrated experimentally by Clanet and Searby jllj . 
The mechanism [TT] describes the early stage burning of a 
flame ignited at the center of the end face of a closed tube, 
leading to the transition of an initially hemi-spherical 
flame kernel into a flnger shaped front as illustrated in 
Fig. [Tj The tip of this finger fiame experiences short but 
powerful acceleration until the flame skirt touches the 
tube walls. Thereafter, the flame acceleration stops, the 
flame skirt catches up with the tip rapidly and the flame 
inverts into a "tulip" shape. This finger flame accelera- 
tion phenomenon was flrst studied in the context of the 
"tulip flame" formation [llj . However, it was pointed 
out in Refs. [51 [H] that the notion of a "tulip flame" 
is too ambiguous, because it may be attributed to any 
concave flame front with a cusp pointing to the burnt re- 
gion, as demonstrated by several combustion phenomena 
of different physical nature [T5H17) . In particular, accel- 
erating turbulent flames in tubes/channels with no-slip 
at the walls, laminar [8l|9] or turbulent [5 , also exhibit an 
irregular "tulip" flame shape. Furthermore, a tulip-like 
flame shape is also relevant to oscillating flames [T5ti23j . 
To avoid any ambiguity, and recognizing that oscillat- 
ing flames imply non-slip at the walls and relatively long 
flame propagation [51 [53] , in the present paper we focus 
only on the the laminar finger fiame acceleration during 
the initial stage of burning, without considering other 
manifestations of tulip fiames. The acceleration appears 
to proceed in a clearly exponential manner, as shown in 
a number of works [TTl [T^ [53]. This is an important 
effect for subsequent DDT, since powerful precursor ac- 
celeration may create a leading shock wave responsible 
for pre-heating of the fuel mixture. 

A quantitative theory of finger flame acceleration in 
cylindrical tubes was developed in Ref. [T5] by assuming 
incompressible flow. The theory shows that the maxi- 
nrum velocity in the laboratory reference frame achieved 
by an accelerating flame tip is (29 — 1)QSl, where Sl 
is the planar unstretched flame speed and Q = Pu/pb 



the initial ratio of the fuel mixture density to the burnt 
gas density. For slow hydrocarbon flames with Sl ~ 40 
cm/s and 8 « 8, this yields a maximum velocity of the 
flame tip ^ 1205*^ ^ 48 m/s, which is still consider- 
ably lower than the sound speed. The situation, how- 
ever, becomes quite different for fast hydrogen-oxygen 
flames, with 5^ ;< 8 m/s, O ~ 8, and the "incompress- 
ible" estimate of Ref. [H] yields (26- l)eS'L ~ 960 m/s, 
which exceeds the sound speed in the mixture, 530 m/s, 
and thereby fundamentally violates the incompressibil- 
ity assumption. Consequently, in order to describe flnger 
acceleration of fast flames properly, gas compressibility 
needs to be accounted for. 

We next note that while the theoretical analysis pre- 
dicts exponential flame acceleration with the incompress- 
ibility assumption [TMlOj. various experiments showed 
moderation of the initial exponential regime with time 
and the possibility for the flame tip velocity to satu- 
rate to a supersonic speed in the laboratory reference 
frame [31 [H E5H57] . Existence of such a saturation ve- 
locity, which correlates with the Chapman- Jouguet de- 
flagration speed, follows from the basic theory of defla- 
gration and detonation fronts [281 ES] ■ Recent numerical 
simulation and analytical theory demonstrated the same 
tendency: gas compressibility moderates the initial expo- 
nential acceleration of the flame it to a slower one, leading 
to saturation of the flame speed [TOl I3TH33] . In addition, 
simulations [7] showed that in obstructed channels fast 
flames with a relatively high initial Mach number exhibit 
noticeably lower acceleration rate as compared to slow 
flames. Since acceleration of finger fiames has much in 
common with ultra-fast flame acceleration in obstructed 
channels [7] , a similar influence of gas compressibility on 
the flnger flame acceleration is expected. 

The purpose of the present paper is to study finger 
fiame acceleration analytically, experimentally and com- 
putationally for various values of the unstretched laminar 
fiame velocity, thus focusing on the influence of gas com- 
pressibility. The early stages of burning in tubes with slip 
adiabatic walls are considered. We first developed an an- 
alytical theory of flame acceleration in 2D (planar) and 
axisymmetric geometries through small Mach number ex- 
pansion up to first-order terms, demonstrating that gas 
compression reduces the acceleration rate and moderates 
the finger flame acceleration noticeably. We then conduct 
experiments for hydrogen-oxygen mixtures with consid- 
erably large initial values of the Mach number, showing 
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FIG. 1: Geometry of an accelerating finger flame. 
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FIG. 2: Flow close to the channel wall. 



the scaled acceleration rate to be much smaller than that 
observed previously for hydrocarbon flames. We also per- 
formed numerical simulations for a wide range of initial 
laminar flame velocities; the results agree well with the 
experiments as well as the theory in the limit of small 
gas compression (small initial flame velocities). Similar 
to previous works, the numerical simulations show that 
the finger fiame acceleration is followed by the formation 
of a "tulip" flame shape, which indicates the end of the 
early acceleration process. 

The paper consists of six sections. In Sections |2 A[|2B[ 
and |2 C| we develop the theory of flame acceleration in the 
early stage of burning. Details of the numerical simula- 
tions are presented in Section[3] Section 4 A describes the 
experimental results. In Section |4B| we compare results 
from the theory, simulation and experiment, followed by 
the conclusions. Resolution and laminar velocity tests 
are presented in the Appendix. 
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FIG. 3: Time limits of the flame acceleration. 



wall of the chamber (^ — 0) from the axis (77 ~ 0) to the 
sidewall (ry = 1), as shown in Fig. [2] The flame separates 
the flow into two regions of fresh mixture and burnt gas. 
Assuming incompressibility for substantially subsonic 
flame propagation, the continuity equation is given by 



dw dv 



0. 



(2.1) 



The boundary conditions are u = at the end wall, C = 0, 
and w = at the side wall, 77 = 1. The flow in the 
fresh mixture (labeled "u") is assumed to be potential, 
so Eq. (lO yields 



2. THEORY OF FINGER FLAME ACCELERA- 
TION 

We consider a flame front propagating in a 
tube/channel of radius/half-width R with an ideal slip 
adiabatic walls as shown in Fig. [l] One end of the 
tube/channel is kept closed, and the embryonic flame is 
ignited at the central point of the closed end wall. It 
was explained in Refs. J1HI12J that in the axi-symmetric, 
cylindrical configuration a flame front develops from a 
hemi-spherical shape at the beginning to a finger shape. 
Here we first present a 2D (planar) counterpart of this 
formulation, assuming flow incompressibility, and then 
we extend the formulation, both in the planar and axi- 
symmetric configurations, to the case of finite, but small 
compressibility. 



A. Finger flame acceleration for planar geometry 

In the theory we employ the standard model of an 
infinitely thin flame front propagating normally with 
the speed Sl, and use the dimensionless coordinates 
ivA) = {Xt'^)/Ri velocities (w,v) = {uxtUz)/ Sl, and 
time T — SLt/R- A 2D fiame, ignited at the point 
(^1 C) = (0:0)5 is initially semicircular, but the flame 
shape changes as the flame-skirt 7]f moves along the end 



Ci^, 



Ci (1 - ry) , 



(2.2) 



where the factor Ci may depend on time, but is indepen- 
dent of the (radial) coordinate. Recognizing that while 
the flow in the burnt gas (label "b" ) is rotational in gen- 
eral because of the curved flame shape, we can never- 
theless treat it as a potential flow close to the end wall, 
where the flame front is locally planar, see Fig. [2] Sub- 
sequently, Eq. (2.1 1 with the boundary condition at the 
channel axis, u" = at 77 = 0, yields the velocity distri- 
bution in the burnt gas in the form 

Vb = C2L Wb = -C2V- (2.3) 

The matching conditions at the flame front, t] ~ rjf, are 



dr 



1, 



Wb 



Vb, 



e-1. 



(2.4) 



(2.5) 



(2.6) 



Here Eq. (2.4) specifies the fixed propagation velocity Sl 



of the flame front with respect to the fuel mixture, while 
Eqs. (2.5) and (2.61 describe continuity of the tangential 



velocity at the front and the jump of the normal velocity. 



respectively. It is noted that the condition (2.5) apphes 



only at the flame skirt close to the wall. Substituting 
Eqs. (|2:2l)-(l231) into Eqs. (pl|-p:6|, we find Ci = C2 = 



(O — 1). Consequently. 



Vu = Vb = {Q - 1)^, 

wu - (e-i)(i-77) 

Wb = - (6 - 1) ?7, 



(2.7) 



and the evolution equation for the flame skirt, Eq. (2.4 1, 
becomes 



dT 



-{Q-l){l-Vf) = l, 



(2., 



which can be integrated with the initial condition rjf — 
at r = as 



Vf 



e 



e-i 



{l-exph(e-l)r]}, 



1 / e-1 

T = —— T In I 1 — — rjf 



e-1 



e 



(2.9) 



According to Eqs. (2.8) and (2.9), we identify two 



regimes of flame propagation, namely, those with the 
flame skirt close to the axis and the wall {rjf ^ 1 and 
1 — rjf <C 1, respectively). In the limit of 77^- <C 1 and 
{Q — 1)t ^ 1, the flame propagates as drif/dr — 9, 
rjf = Or (i.e. Rf — QSl, Rf = BS^i), which is re- 
lated to the expansion of a semicircular flame front. In 
the limit of 1 — 77/ <C 1, a locally planar flame "skirt" 
approaches the wall, the radial velocity of the fresh fuel 
mixture tends to zero, and the flame skirt propagates 
with the planar flame speed with respect to the end wall 
of the channel, drjf /dr — 1 (i.e. Rf = Sl)- The time 
of the transition from the hemi-circular to the "flnger"- 
shape flamefront can be estimated as 



1 



'^sph 



e-1' 



when the position of the flame skirt, Eq. (2.9), is rjf^sph 



'l-e- 



e/ (e - 1) w 0.63 e/ (e - l) sothe transition 
occurs approximately when the flame skirt has moved 
more than half-way to the side wall of the channel. Sub- 
stituting rjf — 1 into Eq. (2.10) we flnd the time when 
the flame skirt touches the tube wall 



'wall 



ine 
e-1' 



(2.11) 



The results (2.10) and (2.11 ) are shown in Fig. Islas func- 
tions of the expansion factor e. It is clearly seen from 
Eqs. ( 2.10[ ) and (2.11) and Fig. [slthat the acceleration 
[Tsph < Twaii) occurs if B > 6. For e = 5 ~ 10 we have 
T,ph « 0.11 - 0.25, while Ty,aii = 0.25 - 0.4. 

To determine the evolution of the flame tip, we con- 
sider the flow along the channel axis, 77 = 0, as shown in 
Fig. |4] The flame shape is considered to be locally planar 
in the vicinity of the centerline (77 — > 0) . In that limit the 
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FIG. 4: Flow close to the channel axis. 



*^ 



flow can be considered to be potential, with the longitu- 
dinal velocity component v determined by an equation 
similar to Eq. (2.3 1. The solution for the longitudinal 



velocity component v has to coincide with Eq. (2.3) at 
7; -^ 0, ^ ^- 0. Thus, the longitudinal velocity compo- 
nent in the burnt gas Vb along the centerline (7; = 0) is 



governed by Eqs. (2.3 1 and (2.1). We stress that this rea- 



soning does not hold away from the axis in the burnt gas, 
where the flow is rotational. Still, only the gas velocity 
along the centerline is utilized in the present analysis. 
Based on the condition of a flxed propagation velocity of 
a planar flame front with respect to the burnt gas. 



di 



tip 



dT 



- ufc = e. 



(2.12) 



with Vb given by Eq. (2.7), we arrive at the differential 



equation for the flame tip, 
d£,tip 



dT 



(e - i)iup = e. 



(2.13) 



with the initial condition (^tip (0) = 0, and the solution 

e 



6^p=Q^{exp[(e-l)r]-l}, 



(2.14) 



which yields a semicircular flame front just after igni- 
tion, (e — 1)t ^ 1, when Eq. (2.14) is reduced to 



(2.10) ^tip = er = 77/, and the transition to the exponential 



acceleration thereafter. According to Eq. (2.14), the 



growth rate during the exponential stage of acceleration 
for the planar geometry is given by 



'^0,pl 



e-1. 



(2.15) 



At the end of the acceleration, when the flame skirt 
touches the wall, we have t — T^aii, so the position of 
the flame tip is 

iwall = itip {Twall) = B, (2.16) 

or Z^aii = ei? in dimensional units. Therefore, account- 



ing for Eq. (2.13), we have 



de 



tip 



dT 



Utip 



= e\ 



(2.17) 



To describe the flame shape and evaluate the total in- 
crease in the flame surface area during the flame accel- 
eration, we assume, realistically, that at the end of the 



acceleration the flame shape is almost self-similar, with 
ry/ ss 1 and ^up « exp [(9 — 1)t]^ 1. We therefore look 
for the flame shape in the form 



= Cupir)- f{T,ri), 
f{T, rj) ^ if (rj) exp [(6 - 1) r] 



(2.18) 



The accuracy of such an approximation is Q^^ ^ 1, 
which is acceptable for typical flames. Then the equa- 
tion of flame evolution with respect to the burnt gas can 
be written as 



5/ 
dr 



■Wb 



5/ 
drj 



Vb^Q 



dl 
drj 



1/2 



e|. (2.19) 



Accounting for the exponential state of flame accelera- 
tion, Eq. (2.181 and the velocity distribution (2.3), we 
reduce Eq. p.im to 



|e-(e-i),i|.e, 



(2.20) 



with the boundary condition at the axis, ip (0) — 0, yield- 
ing the solution 



'fiiv) 



8 



e-1 



In 1 



e-1 



e 



-T] 



(2.21) 



with 



iv, t) 



e 



e-1 



1 - In 1 



e-1 

exp[(e-l)r]-l}. (2.22) 



Then the maximum increase in the flame length, achieved 
when the flame skirt touches the wall, is 



L^/2R= f Jl 



dl 
drj 



drj 



dl 
dr] 



dl] = 



= /(I, T^all) = 95 (1) exp [(e - 1) T^all] = 



e^ine 
e-1 ■ 



(2.23) 



B. Influence of gas compressibility on acceleration 
rate for planar geometry 



As long as Ma ^ 1, we can treat the flow ahead of 
the flame as isentropic. In this case, to first order in 
Ma, we have the following relations for the scaled density 
Pu = Pu/Pu,o, pressure P„ = Pu/Pufl, and temperature 
Tu — Ty^/T^xi in the fuel mixture 



1 + 



7- 1 u„ 



1 



2 cs,o. 
1 + Ma{Q - 1) {^up + 1) : 



csfi 



(2.25) 



P„ 



1 



7- 1 «u 



1-1-7- 



2 csfi J cs,o 

1 + Ma-i{Q - 1) {iup + 1) , 



(2.26) 



Tu = 1 + ^^ «l-h(7-l) « 

V 2 csfij csfl 

« l + Ma(7-l)(e-l)(6»p + l), (2.27) 

where 7 = Cp jCy is the adiabatic index. It is noted that 
Eqs. (2.25)-( 2.27[ ) also define the rigorous mathematical 
limit of validity for the present theory, i.e. v^/csfi ^ 1- 
The matching relations at the fiame front are 



~ I d£,tip \ . I d£,tip 



(2.28) 



Pu + Pu 



~ I d£,tip 



V dr 



P 



bf + Pbf I 



/ d£,tip 



dr 



Vbf 






T,, 



bf 



1 fd^tip 
2Cpno V dT 



Vbf 



(2.29) 



(2.30) 



where Q = Q/CpTufi = e — 1 is the scaled reaction heat 
release. With the first-order approximation for small 
Ma, we reduce Eqs. (12. 291) and p^ to 



Pu = Pbf, Tu + e-l^Tbf. 



(2.31) 



We next consider the problem, with first-order ac- 
curacy for the initial flame propagation Mach number 
Ma = Sl/cs,o ^ 1, where cs,o is the initial sound speed 
of the fresh mixture. This approach is conceptually close 
to that of Refs. [32 ) I33 j . The compressible counterpart 
of Eq. (2.12) for the dynamics of the flame tip becomes 



d£,tip 



"bf 



^, 



(2.24) 



where Vbf — Vb{£,f, t) is the flow velocity of the burnt gas 
in the ^-direction just at the flame front, and i? = Pu/ Pbf 
the instantaneous gas expansion factor, with 'Bq = Q. 



Usin g the perfect gas law, puTu = pbfTbf, as well as 
Eq. (p^, we find 



d = l- 



e-1 



e-Afa(7-l)(e-l)2 (e,,p + 1) . (2.32) 



While, according to the Euler equation 

dvb ^ dvb 1 dPb 

'^T ^'"b^^y = - — ^^^ (2-33) 

OT dS, Pb d^ 

pressure is uniform in the burnt gas, Pbir) = Pbf, up 
to the first-order in Ma <C 1, it however grows in time, 



and thereby increases the temperature and density of the 
burnt gas due to adiabatic compression. 

We next consider propagation of the nearly planar 
flame "skirt". Within the accuracy of Ma, the pressure 
in the fuel mixture between the flame front and the side- 
wall is the same as that in the burnt gas, P, = Pb{T) — 
P(,f = Pu- Thus the density and temperature of the fuel 
mixture around the flame skirt are the same as those of 
the fuel mixture just ahead of the flame tip, since in both 
cases we have adiabatic compression and the same final 
pressure. Consequently, the continuity equation for the 
fuel in the domain between the flame shirt and the side 
wall, r]f < 7] < 1, takes the form 



dws 1 dps 1 dPs 1 dPu 

dri ps dr jp^ dr jp^ dr 



(2.34) 



with the following solution satisfying the matching rela- 
tion at the side wall, 77 = 1, w^ = 0, 



1 dPu ,^ . 
1 dPu ^^ . 



(2.35) 



Here the subscript "s, /" designates the flow velocity just 
ahead of the flame skirt. Using the matching relations at 
the flame front, we find the velocity in the burnt gas at 
the flame skirt, subscripted by "s,6, /", as 

1 dp 

Ws,bj = Wsj-^ + 1 = -^-7^ (1 - ?7/)-i? + l- (2.36) 

We subsequently solve the continuity equation, equiv- 



alent to Eq. (2.34 1, in the burnt gas with Eq. (2.36) to 
find 

1 dPu 1 dPu 



iPu dr 



iPu dr 



1 dP, 

iPu dr 



J^(l-77)-i?-M. 



(2.37) 



Similarly, the continuity equation for the burnt gas 
around the symmetry axis takes the form 



dwb dvb 



1 dp, 



dr] dS, 7P„ dr 

which can be integrated as 



(2.38) 



Vb = i-&^1)^- 



1 dP 



iPu dr 



-i, Wb ^~{d~l)T]. (2.39) 



Consequently, the evolution equation for the flame tip, 
Eq. ( |2.24| , becomes 



'^'- = (^-i-^^U..+^. 



dr 



iPu dr 



(2.40) 



Fi nally substituting Eqs. ( |2.26[ ) and ( |2.32[ ) into 
Eq. (2.401, neglecting the second- and higher-order terms 



in A/a, and accounting for the zeroth-order approxima- 
tion, Eq. (2.14), we obtain 



d£,tip _ , ,, .„ ^s2 ^2 



dr 



Maj (e - 1)^ eL + ^i,pi6^p + ei, (2.4i) 



with 



^upi - (e - 1) [1 - Ma (e + 2 (7 - 1) (6 - 1))], (2.42) 



ei = e-A/a(7-i)(e-i)^ 



(2.43) 



In the limit of incompressible flow, Ma = 0, we have 
<^i.pi = Q — 1, 61 = Q, and Eq. (2.41) fully reproduces 



Eq. (2.13). Accounting for gas compressibility, we obtain 



moderation of the flame acceleration in Eq. ( 2.41 1, which 



is described by two types of terms: linear and nonlinear 
with respect to ^up- The linear term does not change 
the exponential state of the flame acceleration, though 
they reduce the acceleration rate to cri^pi as compared to 
0—1 for the incompressible flow. At the very beginning, 
for £,tip -^ 0, the flame acc elerat ion is moderated by the 
linear terms only, and Eq. (2.41) reduces to 



d£,tip 
~d^ 



0'l,pl£,tip + 01, 



with the solution 



■,tip 



9i 

cri,pl 



[exp ((Tl^piT) - 1]. 



(2.44) 



(2.45) 



The nonlinear term of Eq. (2.41 ), however, becomes im 



portant quite fast and modifles the exponential state of 
flame acceleration to a slower one. The complete analyt- 



ical solution to Eq. (|2.4ip is given by 

t.tip 



201 [exp {a2T) - 1] 



(0-2 - cri,p/) exp (cr2T) + (0-2 + (Jl.pl) ' 



(2.46) 



where a2 ^ Jcrl . + AMa-yOi (0 - 1)' 



C. Influence of gas compressibility for the axi- 
symmetric geometry 

Now we reconsider the problem in an axi-symmetric, 
cylindrical tube, Ref. |12j . incorporating compressibility 
with the accuracy of the first order for fiame propagation 
Mach number. 

The incompressible continuity equation in the axi- 
symmetric geometry takes the form |12j 



ld_ 
77 dr] 



{r]w) 



dv 



(2.47) 



and the axi-symmetric counterparts for the velocity dis- 



Eqs. (2.8h(2.9), are given by 



tributions, Eq. (2.7), and the evolution of the flame skirt, 



f, - 2 (9 - 1) r,/e 
1 



Wb 



(9- 1)77/77, 



*^-(9-l)(l-77P = l 



(2.48) 



(2.49) 



This result exceeds considerably (by a factor of about 
2) its 2D counterpart (9 — 1), see Eq. (2.151, and it is 



slightly smaller than the model estimation 29 of Clanet 
and Searby ill_. 

Now we account for small, but finite gas compression. 
With the axial velocity «„ given by Eq. (2.48), the axi- 
symmetric counterparts of Eqs. ( 2.25 )-( 2.27) and (2.32) 
are 



7- 1 w^ 



P« = 1 



2 Co, 
1 + Ma i2a^up + 9 - 1) , (2.57) 



9 



where 



Vf = 


— tanh (ar) , 
a 


^wall -- 


1 , /e + a 

— In 

2a \9 — a 


Tsph = 


1 
2^' 


a — 


V0(e-i). 



(2.50) 



(2.51) 



It can be readily shown from Eq. (2.50) that ac- 
celeration is possible (i.e. Tgph < T^aii) if 9 > 

1 - [(e - l)/e + l)f)~^ « 1.27. For 9 = 5 - 10, Eq. 

2l^ yields Tsph « 0.05 ~ 0.11, and r^aii ~ 0.19 ~ 0.32. 
These quantities are much smaller than those of the pla- 



nar geometry, Eqs. (2.10) and (2.11 ), making an indirect 



proof that for the axisymmetric geometry acceleration 
proceeds faster th an th e planar one. 

With the result (2.48), equation for the flame tip, £,tip — 
Vb — Q, becomes 



d^tip 

1l^ 



2(9-l)77/(r)6,p(T) = 9, 



d£,tip 
with the solution 



2a tanh (ar) S,tip — 9, 



e 



tip 



2a 



sinh (2aT) . 



(2.52) 



(2.53) 



(2.54) 



which also yields ^waii = 9 similar to the 2D re- 
sult (2.16). At sufficiently late times we have rjf Ri 1 
and Eqs. (|2.52|)-(|2.53|) reduce to 



d^tip 



9 



9 



2a 



(2.55) 



so the fiame tip accelerates almost exponentially, with 
the acceleration rate 



Co,, 



2a = 2^/6(9-1) 



(2.56) 



1 



7- 1 ^'n 



Tu = [l 



2 Co. 
1 + Maj {2a^Hp + 9 - 1) , (2.58) 



7- I v^ 
2 Co. 
1 + Ma (7 - 1) {2a(trp + 6 - 1) , (2.59) 



7? = 9 - Ma{j - 1)(9 - 1)" ( 2-^up + 1 



Following the strategy of Section |2B[ we find 

Id dVb 1 rfP„ 

-^(77Wfc) + ^- = fT-T' =^ 

r]df] d^ 7P„ dr 



Vb 



Wb 



2(79-l)77K-4r^?, 



iPu dr 



{■d - 1) rifij. 



(2.60) 



(2.61) 



(2.62) 



Similar to Eqs. (2.49) and (2.50), the flame skirt position 
is given by 

(2.63) 
where a = ^Jd [d — 1), and the evolution equation for 
the flame tip is 



1 



rjf — — tanh (aT) , 



d^tip 



2(7?- 1)77/ -;^^jW + ^. (2.64) 



Holding the zeroth- and first-order approximations for 
Ma in Eq. ( |2.64 ), we rewrite it in the form 

d£,tip 



dr 



= [2a tanh (ar) - 2Ma-da] £,, 



tip 



-AMaaa tanh (ar) ^. 



tip 



(2.65) 



or 



-|^ = {2ai tanh (air) - 2Ma a [1 + 7 (9 - 1)]} £.tip ~ 

CLT 

( arB 1 

-2Ma i 2 h (2a2 - B) tanh (ar) } ^f^ + 9i 

t cosh (ar) J 

(2.66) 



where 9i is the same as in the planar geometry, see 
Eq. ^A3\, and 



ai = v/ei(ei-i), 
B =. (7-i)(e-i)(2e- 



1). 



(2.67) 



In general, Eq. (2.66) has to be solved computation- 



ally, but we shall integrate it analytically with several 
asymptotic approaches. First, in the limit of Ma = 0, 
we have 0i = O, cui = a and Eq. (2.661 reproduces 



Eq. (2.53). Even for finite gas compressibility, the effect 



of the nonlinear term in Eq. (2.66) is negligible in the 
very beginning, hence Eq. (2.66) can be approximated 
by 

'^^*'P - {2aitanh(air)- 

2Maa[l + 7(e-l)]}6zp + ei, (2.68) 



dT 



with the solution to the first-order approximation for Ma 
being 



6 



ti'p 



2ai 



sinh (2Q;ir) 

e 



2Ma— [H- 7 (6 - 1)1 cosh (ar) In cosh {(xt) 
a 



(2.69) 



The nonlinear term of Eq. (2.66), however, becomes im 



portant quite fast, modifying the state of flame accel- 
eration to a slower one, and hence making the asymp- 
tote (2.69) incorrect. However, at a sufficiently late stage 
of the acceleration, we can approximate tanh (air) ^ 1, 
hence r]f ~ '&/a, and Eq. (2.66) is reduced to a form 



similar to that of the 2D, Eq. (2.41) , 



di, 



tip 



dT 



where 



= -Ma'4>^%p + ai^axi£,tip + ©1, 



V- = 2 (6 - 1) (267 - 7 -H 1) , 



(2.70) 



(2.71) 
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FIG. 5: The sketch of the grid with variable resolution used 
in numerical simulations. 



acceleration rate from (jQ^axi fo '^i.axij the linear terms 
do not change the exponential state of acceleration while 
the nonlinear term of Eq. (2.70) modifies the exponential 



state of flame acceleration to a slower one as soon as it 
becomes important. 



It is emphasized that the result (2.73) does not re- 



produce the asymptote (2.69) as they are related to op- 



posite limiting cases. Consequently, it is expected that 
for a certain range of Ma the complete (numerical) solu- 
tion to Eq. ( 2.66[) lies in between the two values given by 
Eqs. KM and |2.73). 



3. NUMERICAL METHOD, BASIC EQUATIONS, 
BOUNDARY AND INITIAL CONDITIONS 

We perform numerical simulations of the hydrody- 
namic and combustion equations including transport 
processes (thermal conduction, diffusion, viscosity) and 
chemical kinetics in the form of Arrhcnius equation. Both 
2D planar and axisymmetric cylindrical flows are investi- 
gated. In the general tensor form the governing equations 
are given by 



| + ii|:{"'V.) + |(...) = o. ,3.1) 



O"!, 



<^0,i 



1- Ma 



49-1 

e + -.^ (7 - 1) (e - 1) 



26 



and with the solution, Eq. (|2.46|) 



2ei[exp(f72T)-l] 



((72 - (Ji^axi) exp (cr2T) + ((72 + Cri,aa;i) 



(2.72) 
(2.73) 



AQiMa^p in the axisymmetric con- 



where (72 = 

figuration. 

recovers the properties of its planar counterpart: in the 

limit of Ma = we have ai^axi = '^o,axi, and Eqs. (2.70) 



Obviously, the resuh (2.70) to (2.73) fully 



and (2.73) reduce to Eq. (2.55); accounting for gas com- 



pressibility, we obtain linear and nonlinear moderation of 
the flame acceleration with respect to S,tip- Reducing the 



dt 



{pUr 



dt 



dt 



{PUZ 



de Id 



d dP 

{pUrU^ - (r,z) + ^^ + i'P = 0, 



d 



dr 



1 d 



[r'^ {pUrU^ - C,r^z)\ + 

9,2 . \ dP 
dz ' dz 



r^ dr 
d 



(3.2) 



(3.3) 



r'^ ((e + P)Ur- Cr^rUr - Cr.z^z + Qr)] + 



+ -7r[{£ + P)Uz-Cz,zUz-Cr,zUr + qz]^0, (3.4) 

dz 



dt 



(PY) 



1 d 

ri dr 


r-T puiY-——- 
\ Sc or 


oz \ Sc oz 1 



= -^ cxp{- Ea/RT) 
TR 



(3.5) 



where /3 — and 1 for 2D and axisymmetric geometries, 
respectively, 

e^p{QY + CvT) + ^{ul + ul) (3.6) 

is the total energy per unit volume, Y the mass fraction 
of the fuel, Q the energy release from the reaction, and 
Cv the heat capacity at constant volume. The energy 
diffusion vector qi is given by 



qr = -M 



-^J■ 



CpdT Q_dY 

Pr dr Sc dr 

C^dT Q_dY 

Pr dz Sc dz 



(3.7) 



(3.8) 



In the 2D configuration (/? = 0) the stress tensor C^j 
takes the form 



G, 



dui duj 2 duk 
dxi 



dxi 3 dxk *'"' 



(3.9) 



where R « 8.31 J/(mol • K) is the universal gas constant. 
We consider a single-step irreversible reaction of the first 
order with the temperature dependence of the reaction 
rate given by the Arrhenius law with an activation en- 
ergy Ea and the factor of time dimension tr. In our 
simulations we took Ea/RT^ = 32 in order to have bet- 
ter resolution of the reaction zone. The factor tj^ was 
adjusted to obtain a particular value of the planar flame 
velocity Sl by solving the associated eigenvalue problem. 
The flame thickness is defined as 



^/ 



fJ-u 



PtPuSl 



(3.15) 



where pu = 1.16kg/m'^ is the unburned mixture den- 
sity. It is noted that Lf is just a mathematical parame- 
ter of length dimension related to the flame front, while 
the real effective diffusion flame thickness is considerably 
larger [521 US] . We took initial temperature of the fuel 
mixture T„^o = 300ii', initial pressure P^fi — 10^ Pa, 
speciflc heat ratio 7 = 1.4, and 6 = 8. We performed 
the simulations for a rather wide range of initial Mach 
number Ma = Sl/cs^o = 10"^ -- 1.6 x 10"^, with 
the lower and upper values being relevant to hydrocar- 
bon and hydrogen-oxygen flames, respectively 127! . We 
used the tube diameter 2R = 150L/ and channel width 
2R = lOOLf for axisymmetric and 2D simulations, re- 
spectively. 

Similar to the theoretical analysis, we adopt slip and 
adiabatic boundary conditions at the tube walls: 



while in the axisymmetric geometry {f3 — 1) it reads 

(3.10) 



2p ( n'^Ur dUz Ur 

3 \ dr dz r 



^"'"" -i \ dz dr 



Sr,j 



dUr 



duz 
dr 



Finally, the last term in Eq. (3.2 1 takes the form 



2p f Ur duj. duz 
3 V r dr dz 



(3.11) 



(3.12) 



(3.13) 



if /3 = 1, and "0^ = if /? = 0. Here p, is the dynamic vis- 
cosity, and Pr and Sc the Prandtl and Schmidt numbers, 
respectively. 

We take unity Lewis number Le = Sc/Pr — 1, 
with Pr = Sc = 0.75; the dynamical viscosity is 
^ = 1.7 X 10~^Ns/m^. The fuel-air mixture and burnt 
gas are perfect gases with a constant molar mass m — 
2.9 X lO-^kg/mol, with Cy = 5R/2m, Cp = 7R/2m, 
and the equation of state 



P = pRT/m, 



(3.14) 



n-u = 0, n-Vr = 0, 



(3.16) 



where n is the unit normal vector at the walls. At the 
open face end of the tube/channel non-rcflccting bound- 
ary conditions are applied. As initial conditions, we used 
a semi-circular flame "ignited" at the channel axis at the 
closed end of the tube, with its structure given by the 
analytical solution of Zel'dovich and Frank-Kamenetskii 

mm 



T = T,, + in - n) exp (-T^^T^/L/) , 



if z'^ + x^ < r? 



T = er„, if z2 



2 ^ 2 

■ X > rr 



Y^in- T)/{n - Tf), 

Urr. =0, Uz = 0. 



P = Pu, 



(3.17) 



(3.18) 



(3.19) 



Here r/ is the radius of the initial flame ball at the closed 
end of the tube. The finite initial radius of the flame 
ball is equivalent to a time shift, which requires proper 
adjustments when comparing the theory and numerical 
simulations. 

The simulations used a 2D hydrodynamic Navier- 
Stokes code adapted for parallel computations [37]. The 
numerical scheme is second order accurate in time and 




FIG. 6: Schlieren images of finger fiame propagation for stoicfiiometric H2/O2 mixture at pressure 0.2 bar. The shown images 
are evenly distributed in time with interval of 100 ms. Note that the shift of image's left boundary is 10 mm to the right as 
compared to the position of the end wall and ignition point. 



fourth order accurate in space for convective terms, and 
second order in space for diffusive terms. The code is 
robust and accurate; it was successfully used in aero- 
acoustic applications. 2D and axisymmetric simulations 
were conducted. We used mesh with variable resolution 
in order to take into account the growing distances be- 
tween the tube end, the accelerating flame and the pres- 
sure wave, and to resolve both chemical and hydrody- 
namic spatial scales. Typical computation time for one 
simulation required up to 10^ CPU-hours, hence imply- 
ing the need for extensive parallel calculations. 

A rectangular grid with the grid walls parallel to the 
coordinate axes was used. The sketch of the calculation 
mesh used in simulations of flame acceleration from the 
closed tube end is shown in Fig. [5J To perform all the cal- 
culations in a reasonable time, we made the grid spacing 
non-uniform along the z-axis with the zones of fine grid 
around the flame and leading shock fronts. For majority 
of the simulation runs, the grid size in the z-direction 
was 0.25L/ and 0.5i/ in the domains of the flame and 
leading pressure wave, respectively, which allowed reso- 
lution of the flame and waves. Outside the region of fine 
grid the mesh size increased gradually with 2% change 
in size between the neighboring cells. In order to keep 
the flame and pressure waves in the zone of fine grid we 
implemented the periodical mesh reconstruction during 
the simulation run I36| . Third-order splines were used for 
the re-interpolation of the flow variables during periodic 
grid reconstruction to preserve the second order accuracy 
of the numerical scheme. 
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FIG. 7: Scaled flame tip position versus scaled time for differ- 
ent pressures: P — 0.2 bar (triangles), P = 0.6 bar (squares), 
P = 0.75 bar (upside down triangles). Different color for each 
pressure correspond to different experimental runs. 



4. RESULTS AND DISCUSSION 

A. Experimental results 

Experiments were performed in a channel with rect- 
angular cross-section (50 x 50 mm), 6.05 m long with 



10 



3300 
2300 
1300 

^(!0 



-r 

10 



Z, R 



15 



—r 

20 



FIG. 8: Temperature field evolution and "tulip" formation 
for planar geometry, — 14, and Mach number Ala — 0.005 
at different time instants. Time instants are equally spaced 
in the range of (0.04022 ~ 0A2161)SLt/R. 



24 transparent ports for photo-gauges. A high-speed 
schlieren system with stroboscopic pulse generator and 
high speed camera, germanium photodiodes and piezo- 
electric transducers were used to record the flame evolu- 
tion. The experimental facilities, including ignition con- 
ditions, are described in detail in Ref. [Ml |^ and the 
references given therein. 

The experiments were conducted with the fast- 
burning, stoichiometric hydrogen-oxygen mixture at ini- 
tial pressures of 0.2 ~ 0.75 bar. The mixture was pre- 
pared by precise partial pressure method with deviation 
less than ±0.2% with respect to the H2 fraction. By 
changing the initial pressure we vary the value of Sl 
and the initial flame Mach number. The initial ambi- 
ent gas mixture temperature was Tg — 293K; and the 
corresponding sound speed is 531 m/s. The initial den- 
sity ratio is in the range 9 = 8.0 ~ 8.317, depending 
on the initial pressure. It is noted that for such highly 
reactive fuel mixtures the influence of the ignition source 
on the initial flame dynamics is negligible, as compared 
with the influence of mixture reactivity, expansion ratio 
and tube diameter. The ignition energy of 2-10 mJ is 
about 10'^ times lower than the released combustion en- 
ergy at the initial stage of finger flame propagation (with 
the diameter of flame ball smaller than 1 cm) . 

Figure [6] shows the schlieren images of "flnger" flame 
propagation at a pressure of 0.2 bar. It is seen that the 
flame tip accelerates exponentially at the initial stage 
of flame propagation, in agreement with previous experi- 
ments [11 and theory [12]. For the stoichiometric H2/O2 
fuel mixture used in the present experiment the influ- 
ences of diffusional-thermal cellular and pulsating insta- 



bilities [5] are ruled out. The influence of hydrodynamic 
(Darrieus-Landau) instability is negligible as well due to 
the strong flame curvature observed in the experiment, 
leading to the Zeldovich-type stabilization of the flame 
front perturbations ,34] . Moreover, even if the Darrieus- 
Landau instability is developed, its scaled exponential ac- 
celeration rate would be about unity T, T , which is up to 
an order of magnitude smaller than the finger fiame accel- 
eration rate of Eqs. (|2.15|), (|2.42|, (|2.56|, (|2.72|, thereby 



diminishing any possible role of the DL instability as 
compared to the finger-flame acceleration effects. Con- 
sequently, in the present experiment the flame front ac- 
celeration could be attributed purely to the flnger-flame 
mechanism of exponential acceleration [11] [12] . 

The experimentally obtained evolution of the scaled 
flame tip position for various pressures is shown in Fig. [Tj 
The definition of the laminar burned velocity Sl,c used 
in the present scaling is described further. When the 
flame skirt approaches the wall, the flame acceleration 
weakens, and the flame shape undergoes transition from 
a convex "flnger" to a concave "tulip" shape. The in- 
stant when the flame skirt approaches the channel wall 
is clearly correlated with the instant at which the expo- 
nential flame tip acceleration terminates, which can be 
seen at Fig. [7j 

Physical parameters of the experiments, as well as the 
scaled acceleration rates obtained from fitting of experi- 
mental results of Fig.lTJfor different pressures, and the ini- 
tial Mach numbers are given in Table [ij The unstretched 
laminar fiame velocity Sl and the thermal ffame thick- 
ness Lx are obtained from the numerical simulation 
of one-dimensional premixed fiame structure employing 
PREMIX code of the CHEMKIN family [HS] with the use 
of updated chemical kinetics mechanism for hydrogen ox- 
idation [3^. The thermal flame thickness Lt is conven- 
tionally defined as [2] {Tb-T^)/ uieix{\dT/dx\), where Tt 
and T„ are the temperatures of burnt and unburnt gases, 
respectively, and max(|9T/9x|) is the maximum of the 
temperature gradient. It is noted that in our case of low- 
pressure stoichiometric hydrogen-oxygen finger flames a 
noticeable pure curvature effect [2] is observed, therefore 
we introduce relevant modiflcations to Sl and Ma, de- 
noted as Sl.c and Mac, in such a mannner that the cor- 
rection parameter C is deflned as C = Sl.c/Sl^ so that 
the estimation for the modifled growth rate is given by 
'^i,C = fi/C. It is seen that, the lower is the pressure, the 
larger is the thermal flame thickness, and consequently 
the more pronounced is the pure curvature effect. 

The pure curvature correction parameter C and the 
modifled growth rate ai^c are estimated as follows. 
Based on experimental observation, we assume, realis- 
tically, that the curvature of the flame tip remains al- 
most constant during the time interval r^ph < t < Ty^jaii, 
with Tsph and Twaii given by Eqs. (2.101 and (2.11), re- 



spectively. Within this time the flame tip radius can be 
approximated as 

Rup = RVf.sph = 0.63 EG/ (6 - 1) « 0.72R, (4.1) 
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P, bar 


e 


Lt, mm 


C 


Sl, m/s 


Sl,c, m/s 


Ma 


Mac 


CTl 


o"i,c 


0.20 


8.00 


2.0 


0.778 


6.855 


5.33 


0.0128 


0.0100 


4.18 


5.373 


0.60 


8.260 


0.5 


0.944 


8.395 


7.928 


0.0157 


0.0149 


4.83 


5.138 


0.75 


8.317 


0.379 


0.958 


8.700 


8.335 


0.0162 


0.0157 


4.90 


5.110 


0.75 


8.317 


0.379 


0.958 


8.700 


8.335 


0.0162 


0.0157 


5.43 


5.660 



TABLE I: Experiment parameters and ai^c obtained by fitting of the results of Fig. JTlfor different pressures. 




FIG. 9: Flame shape evolution and "tulip" formation for axisymmetric geometry, O = 14, and Mach number Ala = 0.005 at 
different time instants. Time instants are equally spaced in the range of (0.02692 -^ 0.26576)5*1,4/7?. Isosurfaces are shown for 
T = 1400^". 



where R is the channel half-width. Consequently, with 



the curvature term in the form V • n 



2/R; 



tip: 



the 



curvature-modified laminar burning velocity at the flame 
tip can be estimated as [2j 



2Li 



Sl,c = CSl = 5l (1 - LtV ■n) = SL[l 



R 



tip 



(4.2) 

As Table |T] shows, modified values of the growth rate 
CTi^c decrease with increasing Ma, while the uncorrected 
growth rate ci increases instead. This demonstrates a 
non-negligible effect of pure curvature for low-pressure 
H2/O2 finger flames, disregarding which could reverse 
the main trend. For this reason, an accurate determi- 
nation of the laminar burning velocity is crucial for the 
analysis of experimental results. This is particularly rel- 
evant for the extraction of the growth rate cti for various 
Ma, since even a 15-20% difference in the laminar burn- 
ing velocity could lead to a completely erroneous conclu- 
sion. 



B. Numerical results and discussion 

We performed numerical simulation of the flame accel- 
eration in tubes with smooth slip adiabatic walls at dif- 
ferent flow parameters. Particularly, we employed planar 
and axisymmetric geometries, and a range of initial flame 
propagation Mach numbers Ma = 10"'^ ~ 1.6 x 10~^. 
We also used two values of the thermal expansion ratio 
6 = 8, 14. 

We illustrate evolution of the temperature field for the 



planar geometry in Fig. [Sjfor 8 = 14 and the initial Mach 
number Ma = 0.005. The simulation plots also show the 
instant of transition from a convex "finger" to a concave 
"tulip" flame shape. It is seen that, similar to the experi- 
mental Figure [6] the flame tip curvature does not change 
significantly from the instant of the transition to the fin- 
ger configuration and until the flame skirt touches the 
wall, thus justifying the assumption of almost constant 
tip curvature made in Sec. |4 A| It should be noted that 
formation of finger-shaped laminar flame fronts in planar 
geometry due to essentially different Schelkin mechanism 
has been also obtained in simulations of premixed flames 
in channels with non-slip walls [311 ESI HO] . Furthermore 
similar shapes of the finger front were observed within the 
context of electrochemical doping in organic semiconduc- 
tors [4TJ |42] , with the electric field playing conceptually 
the same role as the field of gas velocity in the present 
combustion problem. Another interesting physical exam- 
ple of front acceleration and DDT has been encountered 
recently in the studies of spin-avalanches in crystals of 
nanomagnets [121 HI] . 

An axisymmetric counterpart of Fig. [8] is shown in 
Fig. [9| with all parameters being the same except for 
the geometry. Isosurfaces are shown for T = lAOOK. 
The axisymmetic simulation also demonstrates the tran- 
sition from a convex, "finger" flame front, to a concave 
"tulip" fiame, accompanied by a significant reduction of 
the flame surface area and propagation velocity. Similar 
to Figs. [6] and [81 the flame tip curvature remains almost 
constant in a range of scaled times, Tsph < t < T^aiii cor- 
responding to the later stage of finger fiame acceleration. 



In Figure 10 we present the scaled flame tip posi- 
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FIG. 10: Scaled flame tip position versus time for planar geometry, Q 
0.001, 0.0035, 0.005, 0.0075, 0.01: (a) linear scale, (b) logarithmic scale. 



= 14, and Mach numbers Ma = 




200. 



FIG. 11: Scaled velocity of flame tip versus time for 
planar geometry, O = 14, and Mach numbers Ma = 
0.001, 0.0035, 0.005, 0.0075, 0.01. 



tion versus time for the planar geometry with G = 14 
and Ma = 0.001 ^ 0.01 in both linear and logarith- 
mic scales. We see that the "exponential" nature of 
acceleration in the planar case is not very pronounced 
in the early stages, due to the relatively high value of 
Ts-ph ~ 1/(0 ~ 1) compared to its axisymmetric coun- 
terpart Tsj,h ~ 1/(2^/6(6- 1)) [13; the latter being 
approximately twice smaller for high values of than 
in the planar case. Figure [TT] shows the scaled veloc- 
ity evolution for = 14 and the planar geometry for a 
set of Mach numbers in the range Ma = 0.001 ^ 0.01. 
Here we see a significant dependence of the maximum 
flame tip velocity on the initial Mach number. At the 
same time, Fig. [TT] shows that the scaled time of attain- 
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FIG. 12: Scaled maximum flame tip velocity versus Mach 
number obtained in simulations for planar geometry, O = 14, 



is compared to theoretical estimates. Solid line - Eq. (4.8 1 



dashed line - Eq. (4.7 1 



ing the velocity maximum is almost independent of the 
Mach number, which can be demonstrated analytically, 
as follows. If we write the flrst-order correction to S,tip in 
the form ^ up = Ctip.o + Ma ^up.i, where ^tip,o is given by 
Eq. ( |2.14[ ), then Eq. ( |2.32[ ) becomes 

I? = e-Ma{^~l){e-lf{^upfl + Ma^up,i + l)~ 
we- Afa(7 - l)(e - 1) [8 exp ((9 - l)r) - 1] . 

(4.3) 



Similarly, replacing r in Eq. (4.3) by Ti^aii in the form 

Twall = Tyjall.O + MaTy,aU.l, WC find 



e-Ma(7-i)(e-i)(e2-i). 



(4.4) 
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FIG. 13: Scaled flame tip position versus time for axisymmetric geometry, O — 14, and Mach numbers Ada 
0.001, 0.0035, 0.005, 0.0075, 0.01: (a) linear scale, (b) logarithmic scale. 
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FIG. 14: Scaled velocity of flame tip versus time for ax- 
isymmetric geometry, O = 14, and Mach numbers Ma = 
0.001, 0.0035, 0.005, 0.0075, 0.01. 
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FIG. 15: Scaled maximum flame tip velocity versus Mach 
number obtained in simulations for axisymmetric geometry, 
B = 14, compared to theoretical estrimates of Eqs. (4.101 and 
(4ll|. 



Substituting Eq. (4.4) into Eq. (2.11), we obtain the es- 
timation of T-uiaii in the compressible case for the planar 
geometry: 



'wall 



In^ 

ine 
e-1 



l-t-Ma(7-l)(e^-l) 1- 



e-1 
eine 



(4.5) 



For t ypic al « 5 ~ 10 and 7 sa 1.4, the last term 
in Eq. (4.5) can be approximated as ^ Q.2MaQ^. Con- 
sequently, for small Mach numbers. Ma ^ 5/8^ ~ 0.1, 
the quantity Twaii = Ini9/(z9 — 1) only slightly depends 



on Ma, which is substantiated by the simulation results 



of Fig. 11 The weak dependence of T^aii on Ala is con- 
venient for evaluating the total time of the finger flame 
acceleration. We see from the simulation results of Fig- 
ure [11] that the instant of the maximum flame tip ve- 
locity T^aii only slightly increases with increasing Ma, 
remaining in the range T^^aii — 2.2 ^ 2.3. Using the 
simplified analytical expression for r^aii we can also es- 



timate the maximum fiame tip velocity from Eq. (2.41) 



and Eq. (2.44) 



The peculiar feature of concurrent maximums of flame 
tip velocity for various initial Mach numbers can act as a 
useful criterion for the validity of Sl and Ma determined 



14 



in the experiments. While in the numerical simulation we 
set Ma as a parameter and obtain the concurrence of the 
maximum tip velocities mentioned above implicitly, anal- 
ysis of the experimental data can encounter significant 
difficulties due to the uncertainty of Sl, as discussed in 
|4 A[ Since Sl is used in scaling time and velocity, even a 
small inaccuracy in its determination could lead to con- 
siderable shifting of the maximums of the tip velocity 
relative to each other. Thus, if the experimental condi- 
tions imply R « Rup ^ L/, i.e. pure curvature effect 
on flame tip velocity is insignificant at the later stages 
of the flame acceleration, the concurrence of the scaled 
flame tip velocity peaks could serve as an indication of 
correctly determined Sl values for different Ma. 

It is further noted that, similar to the experiments, in 
the present numerical simulations the initial flame veloc- 
ity of a hemispherical flame front is considerably affected 
by the pure curvature effect, which can be seen from Fig. 
11 for the flame tip velocity taken at ISl/R = 0. The 
initial flame radius in the numerical simulation is equal 
4.0L/; thus, similar to Eq. 
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4.2 



we can esti- 



t'O ^init 

mate the correction to the initial flame velocity 5'^^ i^n m 
the planar case, with the curvature term V • n — 1/ Rinit 
in the form: 



Ma 



FIG. 16: Scaled acceleration rate versus Mach number for 
planar geometry, = 8, 14. Values of cti obtained from nu- 
merical simulations are shown by circles for = 8 and by 
squares for = 14. Solid lines correspond to the theoretical 
dependencies give by Eq. (2.421, dashed line - to incompress- 
ible theory, see Eqs. (2.13l-(2.15l, dotted lines - to linear fit 
of simulation data. 



Sl, 



Sl 



1 



L 



f 



Ri 



3 

4' 



(4.6) 



For 6 = 14, the initial tip velocity in the laboratory ref- 
erence frame is Uup/Sl = &SL,init/SL = 0.759 = 10.5, 
which is close to that observed in Fig. [TT] However, the 
effect of pure curvature does not affect the present nu- 
merical simulations considerably, since we have taken rel- 
atively large channel widths, corresponding to i? = 50 Lf 
in the planar case, as specified in Sec. |3J which renders 
pure curvature effects to be negligible for all simulation 
runs. 



for the linear term only, Eq. (4.8), differ significantly. 



with the results of numerical simulation lying closer to 
the exact solution of Eq. (4.7 1. Most importantly, we 



Figure 12 shows the scaled maximum flame tip velocity 
versus the Mach number for the planar geometry and 
G = 14. Analytical estimates for the maximum flame 
tip velocity shown in Fig. 12 are calculated as follows. 



The first theoretical estimate for the maximum flame tip 
velocity accounts for all nonlinear t erms of Eq 
and employs the exact solution, Eq. (2.46), as 



(2.41) 



see that both the simulation and theoretical results show 
signiflcant reduction of the scaled maximum flame tip ve- 
locity with increasing initial Mach number, despite the 
fact that the non-scaled maximum tip velocity could be 
still increasing. Figure [12] shows that the influence of 
gas compressibility is noticeable for the estimation of the 
maximum flame velocity for finger-type flame accelera- 
tion, while previous theoretical studies [11] [12] did not 
account for the significant reduction of the maximum 
scaled tip velocity for high initial Mach numbers. This is 
because they were conducted for slow methane-air flames, 
for which the effect of compressibility did not manifest 
itself. 

We next investigate acceleration of the finger-shaped 
flames at various values of the initial Mach number for 
the axisymmetric geometry. Figure [13] shows the scaled 
fiame tip position versus time for the axisymmetric geom- 
etry for 9 = 14 and the Mach numbers Ma = 0.001 -- 
0.01, in both linear and logarithmic scales. It is seen 
that the stage of exponential flame acceleration is more 
distinctive in the axisymmetric geometry as compared 
to the planar case of FigJTOJ Figure 14 is the axisym- 
metric analogue of Fig. jllj Similar to the planar case, 
the maximum flame tip velocity strongly depends on the 
initial Mach number, but the time T^aiu when the maxi- 
mum flame tip velocity is achieved, only depends on Ma 
slightly. For 9 = 14 used in the s imulation of Fig. [T4| 
we have a « 13.5, and Eq. (2.50) yields r^aii ~ 0.15. 
alytical solution, Eq. (4.7), and the solution accounting The numerical simulation of Fig. |14| shows somewhat 



Ut 



tip, max 

Sl 



-Ma-i{Q~lfeup, 

0'l,pl£,tip,waU + 9i, 



(4.7) 



where S.tip. wall is the flame tip position at time Tj^^aii cal- 
culated from Eq. (2.46), with t = Ti^aii- The second 



estimate is somewhat simplified accounting for the linear 
term only, derived from Eq. (2.44) as 



Ui 



tip, max 

Sl 



= (^i,pi^. 



tip, wall 



+ 9i. 



(4.8) 



"ivith ^up, wall = (9i/cti^p;) [exp(cri,p;T^aH) - !]• It is 
seen from Fig. [12] that predictions of the complete an- 
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FIG. 17: Scaled acceleration rate versus Mach number for 
axisymmetric geometry, = 8, 14. Values of cti obtained 
from numerical simulations are shown by circles for = 8 
and by squares for = 14. Solid lines correspond to the 
theoretical dependencies given by Eq. (2.721, dashed line - to 
incompressible theory (see Eq. ( 2.56 1) , dotted lines - to linear 
fit of simulation data. 



delayed maxiniums of flame tip velocity with t^.q^ ~ 0.2, 
as compared to the theoretical prediction. This delay is 
attributed to the fact that for the axisymmetric case the 
effect of pure curvature is more pronounced in the initial 
stage of finger flame propagation, since the axisymmetric 
counterpart of Eq. (4.6) yields 



Sl, 



1 



2L 



f 



(4.9) 

JL Rinit 2 

However, in the axisymmetric case we have a wide do- 
main with R = 75 L/, see Sec. l3j which renders pure 
curvature effect at the later stage of finger flame acceler- 
ation negligible for all numerical runs in the axisymmet- 
ric case, resulting in concurrent peaks of the flame tip 
velocity in Fig. [T4j similar to the planar case. 

An axisymmetric counterpart of Fig. 12 is presented 
in Fig. 15 Contrary to the planar case, the theoretical 
estimate of the maximum tip velocity shown in Fig. [15] 
is related to two opposite limiting cases descrbed in Sec. 
2C The first estimate is obtained from Eq. (2.68) as 



Ui 



tip, max _ 

Sl 
-2Ma a [1 



{2Q;itanh(aiT^,Qi;) - 
+ 7(e-l)]}6.p,wi + ei, (4.10) 

where the fiame ti p pos ition ^tip,waU at r = r^aii is cal- 
culate d fro m Eq. ( 2.69 ) while ai^axi and Qi are given by 
Eqs. (2.72) and (2.43), respectively. For the second lim- 



iting case of late stage acceleration, the flame tip velocity 
estimate is given by Eq. ( |2.70 ) as 



U< 



tip. max 



-2MaiQ - l)(2e7-7 + lK£p,.„a« + 

+ 0-l,axi£,tip.,waU + Ql, (4-11) 






e 


Sl, cm/s 


Cs m/s 


Ma 


1/i, s-^ 


o"i 


1.0 


8.02 


41.8 


340 


0.00123 


132.6 
129.7 


15.8 
15.5 



TABLE II: Experimental values of cti obtained from Ref. [llj . 
for O ~ 8, low-Mach propane-air mixture and axisymmetric 
geometry. (j> is mixture equivalence ratio, Cs is sound speed 
in unburnt gas mixture, 1/t is experimentally measured tulip 
flame growth rate (for two experimental runs). 



where the flame ti p pos ition S,tip,waU sX t ~ r^aii is cal- 
culated from Eq. (2.73). We see that results of the nu- 



merical simulation in most cases are located in between 



the values given by Eqs. (4.10) and (4.11), as suggested 
in Sec. 2C although, this tendency appears to change 



for higher values of Ma. 



Figure [16] shows the scaled acceleration rate versus 
Mach number for the planar geometry, with = 8, 14. 
It is seen that the agreement between theory and simula- 
tions is significantly better for 9 = 8, indicating that the 
role of gas compressibility increases both with Ma and Q. 
Furthermore, Fig. [16] shows that the growth rate a sig- 
nificantly decreases with increasing Ma, and the isobaric 
theoretical model of the finger-flame acceleration, Refs. 
[TTl 112] , noticeably overestimates the acceleration rate a 
for high-Ma flames. We next note that an accurate es- 
timate of the acceleration rate is of crucial importance, 
e.g., for the analysis of flame-generated shocks, preheat- 
ing of the fuel mixture, pre-detonation run-up distance 
and the DDT onset [301 I3S| • Figure [T7[ shows the scaled 



acceleration rate versus the Mach number for the ax- 
isymmetric geometry, with = 8 and 14. Similar to 
the planar case of Fig. [T6j the agreement of theory and 
simulations is better for 8 = 8 than for Q = 14, which 
indicates a more important role of the compressibility ef- 
fects for larger 8. For 8 = 14 we have even more signifl- 
cant deviation of the numerical data from the theoretical 
predictions as compared to the planar case. Still, the 
trend remains the same: a decreases quite rapidly with 
increasing Ma. 

Figure [18] shows the comparison of the scaled accel- 
eration rates obtained in simulations for the planar and 
axisymmetric geometries, for 6 = 8 in both cases; the 
experimental results for cri obtained in the present work 
and those of Ref. [TT] are also shown. The latter ex- 
perimental data is summarized in Table [ll] for two ex- 
perimental runs for neutrally stable {(f> = 1.0) low-Ma, 
propane-air mixture with 8 = 8.02, and axisymmetric 
geometry of the tube. We see from Fig. [18] that the 
scaled acceleration rate cti obtained in the present exper- 
iments for the high-Afa H2/O2 stoichiometric mixture 
in a channel with a quadratic cross-section are consid- 
erably closer to the planar case than the axisymmetric 
one. Thus, although the flame tip has the hemispherical 
shape, the global flame front acceleration is governed by 
the quadratic cross section of the channel and can be de- 
scribed by the theory of Section 2 A At the same time. 
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FIG. 18: Comparison of scaled acceleration rates versus Mach 
number for experiments and simulations for planar and ax- 
isymmetric geometries, = 8. Values of ai obtained from 
numerical simulations are shown by circles for axisymmetric 
geometry and by squares for planar geometry. Triangles show 
experimental results of Table [l] diamonds - of Ref. iTT], see 
Table III] Solid lines correspond to the theoretical dependen- 
cies given by Eq. (2.421 (planar) and Eq. (2.721 (axisy mmet- 
ric), dashed line - to incompressible theory, see Eq. (2.151 
(planar) and Eq. (2.561 (axisymmetric), dotted lines - to lin- 



end of the finger flame evolution the additional displace- 
ment of the curved elements of the flame would be of the 
order of SLt^aii = SLT^jaiiR/SL ~ 0.25 (in nondimen- 
sional units), i.e. considerably smaller than the flame 
tip displacement ^up = ©■ Due to the elongated flame 
shape resulting from the exponential flame tip accelera- 
tion, the motion of the curved parts of the flame front 
does not considerably influence the flame acceleration. 
Consequently, pressure field variation has only minimal 
influence on the present theory. It is also noted that 
since our problem is evolutionary (non-steady), the max- 
imum flame tip velocity could be more than two orders 
of magnitude higher than Sl, hence implying that the 
low-Mach limit is not applicable, and the pressure field 
cannot be considered steady as, for example, in Bun- 
sen flame problem [i46] . This is beyond the scope of the 
present investigation. 

Finally, it is also noted that the influence of gas com- 
pressibility on the flame acceleration obtained in the 
present paper is qualitatively different from that of the 
DL instability in compressible gases and plasmas [T7I - 
H^ . The compressibility effect renders the DL instability 
much stronger in both linear and nonlinear stages, while 
in the case of finger flame acceleration compressibility 
moderates fiame acceleration in tubes considerably. 



ear fit of simulation data. 



values of cti of Table ITT] [11 for low-Af a, propane-air mix- 
tures are close to the incompressible value ctq for the ax- 
isymmetric case, in agreement with the present theory. 
It is known that flame acceleration in DDT sensitively 
depends on flame confinement and the channel geometry 
[45| . The present numerical simulation and experimental 
results indicate that the geometry of the tube affects the 
growth rate of finger fiamc acceleration significantly, with 
the axisymmetric tube being favorable for faster flame 
acceleration due to the finger fiame mechanism. 

We conclude the present extensive investigation by 
noting that pressure jump at the fiame front could po- 
tentially induce a gas velocity of 0{Sl), and as such 
the pressure field could play an essential role in the 
present phenomena of interest [U US]. Here we recog- 
nize that the present theoretical results are obtained by 
considering only two (quasi-planar) parts of the fiame: 
1) flame skirt near the closed tube end, and 2) flame 
tip in the vicinity of channel axis. Due to the expo- 
nential acceleration of the flame tip, the entire process 
of finger fiame propagation occurs rapidly. For exam- 
ple, for = 8 the characteristic non-dimensional time of 
reaching the maximum flame tip velocity in planar ge- 
ometry is Twaii ~ 0.25. By that instant the flame tip 
has shifted to the non-dimensional position ^up = 6. If 
we assume that the pressure field induces fiame veloc- 
ity of the order of 5*^ at the curved parts of the flame 
(i.e. not at the skirt and the flame tip), then at the 



CONCLUSIONS 



The theory, experiments and numerical simulations of 
the present work show that the growth rate a of the 
finger-flame acceleration from the closed end of a chan- 
nel/tube decreases significantly with increasing initial 
Mach number, Ma. Hence, previous theoretical esti- 
mates of [TI] [12], derived with the incompressible ap- 
proximation, overestimate a for flames with high laminar 
burning velocities, such as the H2/O2 or acetylene/air 
flames. In the present study, we account for gas com- 
pression through expansion for small Mach number up 
to first-order terms, and validate the theoretical analysis 
by numerical simulations and experiment. The present 
numerical simulation and theory show that the maximum 
fiame tip velocity significantly depends on Ma, with the 
scaled time of the maximum flame tip velocity being al- 
most independent of it. The present results collectively 
demonstrate that the geometry of the channel affects the 
growth rate of the finger flame significantly, with the ax- 
isymmetric channel being more conducive for fast initial 
fiame acceleration from the closed end in channels with 
smooth walls. It is emphasized that compressibility ef- 
fects should be taken into account when estimating the 
strength of shock waves generated by the initial finger- 
type flame acceleration, pre-heating of the unburnt fuel 
mixture ahead of the flame, the DDT onset time and 
position. 
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A^//L; 


Umax/ Sl 


/^UmaJSL 


tmaxSL/R 


^tmaxSh/ R 


Zup,*/R 


^Zup.,/R 


1.0 


31.14 




0.3265 




4.5095 




0.5 


34.68 


3.54 


0.3032 


0.2706 


5.282 


0.7725 


0.25 


36.171 


1.43 


0.2973 


0.0059 


5.441 


0.159 


0.125 


36.85 


0.679 


0.2965 


0.0008 


5.474 


0.033 



TABLE III: Resolution tests for planar geometry, Ma — 0.005, O — 
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FIG. 19: Flame tip position versus time for different values 
of the mesh size. Thick solid line correspond to Azf/Lf = 
0.125, dashed line - to Azf/Lf — 0.25, thin solid line - to 
Azf/Lf = 0.5, dotted line - to Azf/Lf = 1.0 



FIG. 20: Velocities of the flame tip versus time for differ- 
ent values of the mesh size. Thick solid line correspond to 
Azf/Lf = 0.125, dashed line - to Azf/Lf = 0.25, thin solid 
line - to Azf/Lf — 0.5, dotted line - to Azf/Lf — 1.0 
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7. APPENDIX: RESOLUTION AND LAMINAR 
VELOCITY TESTS 

In order to check if the adopted resolution is sufficient 
to study the flame acceleration process, we performed 
the resolution tests for the primary results for Ma = 
0.005. The grid size in the flame domain varied between 
0.125L/, 0.25L/, 0.5L/ and IL/. We checked the velocity 
of the flame tip at the instants corresponding to the state 
of the maximum flame tip velocity, as well as the flame 
tip position at the time instant ISl/R = 0.303. The 



resolution test results are presented in Table [111] and in 
Figs. [19] and [20l 



Notation: Azf/Lf is the spatial step in the flame 
grid domain; Umax/ Sl maximum flame tip velocity (see 
Fig. (20l); tmaxSh/R scaled time moment correspond- 
ing to the maximum of flame tip velocity; Zup^^,/ R flame 
tip position at time I^Sl/R = 0.303 (see Fig. (20)). 
^Umax/ Sl increment of Umax/ Sl calculated in the table 
row i as AUmaxii) = Umax{i) - Umax{i - !)■ Increments 
for tmaxSL/R and Zti^^^/R are calculated in a similar 
manner. Resolution in the wave grid domain is equal to 
Azw = 2 X Azf for each run. 

Table [III] and Figs. ]I9] and )20] show good convergence 
of the numerical solution with the increase of mesh reso- 
lution. Resolution tests also showed convergence of time 
corresponding to the maximum flame tip velocity with 
increasing resolution. 

In addition to the resolution tests, laminar flame ve- 
locity tests were performed for Ma — 0.001, 6 = 8, 14. 
The numerical setup was similar to the main numeri- 
cal experiments, except for planar initial flame front and 
tube width 2R = 4L/. That value was chosen to be be- 
low the critical diameter needed for the growth of the 
hydrodynamic instability [34] , so that the flame front re- 
mains planar during the test simulation. For Q = 14, 
the measured laminar planar flame velocity in labora- 
tory frame was Uf lab/SL = 13.97 ± 0.03; for 6 = 8 
Uf, lab/SL = 7.97 ± 0.05. 
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